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Abstract 



Using an ansatz motivated by the classical form of e lr ^, where <f> is the angle 
variable, we construct operators which satisfy the commutation relations of 
the creation-annihilation operators for the anharmonic oscillator. The matrix 
elements of these operators can be expressed in terms of entire functions in 
the position complex plane. These functions provide solutions of the Ricatti 
equation associated with the time-independent Schrodinger equation. We re- 
late the normalizability of the eigenstates to the global properties of the flows 
of this equation. These exact results yield approximations which complement 
the WKB approximation and allow an arbitrarily precise determination of the 
energy levels. We give numerical results for the first 10 levels with 30 digits. 
We address the question of the quantum integrability of the system. 
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In recent years, many interesting questions regarding the quantum behavior of systems 
which classically exhibit sensitive dependence on the initial conditions have been raised, fl 
Classically integrable systems with a perturbation |2| provide rich sequences of transitions 
to chaos when the parameter controlling the perturbation is increased. Recent studies Q 
of the energy spectra of the quantum version of such systems have provided unexpected 
results. Ultimately, the questions raised in this context may become relevant for quantum 
field theory on a space lattice, i.e., a system of coupled anharmonic oscillators. It would be 
interesting to understand if the tool which is used classically to control the effects of the 
perturbation, namely the analysis of small denominators, can also be used for the quantum 
problem. The prerequisite for such a discussion is to have at hand a quantum version of the 
action-angle variables for classically non-linear problems. Due to the ordering problem this 
is a non-trivial issue. 

In this letter, we construct a quantum operator corresponding to the classical quantity 
e 1 ^, cj) being the angle variable, for one of the simplest non-linear problems with one degree 
of freedom, namely the anharmonic oscillator with an energy operator 

H=^ + \muj 2 x 2 + Ax 4 . (1) 
2m 2 

We assume A > 0, but the sign of u 2 is not crucial for the calculations which follow. 
The method proposed can indeed be extended straighforwardly to the case where V is any 
even polynomial bounded from below. Classically, the angle variable satisfies the Poisson 
bracket relation {H, e 1 ^} = —ie l< ^uj(H), where u is one over the derivative of the action with 
respect to the energy. A quantum version of this equation reads 

[H, A] = HAA{H) , (2) 

where A is a function which remains to be specified. In order to construct an operator A 
satisfying Eq. (2), we use an ansatz which in the classical case provides an explicit expression 
for e %< ^ in terms of p and x. The form of the classical ansatz can be obtained directly from 
the well-known expression of the angle variable in terms of the position and the energy. In 
addition, we pick a special type of ordering, namely 

oo 

i = ^ 2 »rt + a„(ff)) . (3) 

n=0 

With this choice of ordering, a solution A multiplied on the right by an arbitrary function 
of H is another solution. In the classical theory, this ambiguity is raised by imposing that 
e 1 ® be on the unit circle. In the case of the potential given in Eq. (1), we obtain a formal 
solution of Eq. (2) provided that for n > 

h 2 - 

—in + 2)(n + l)(2n + 3)K n+2 (H) = [-(2H + hA(H))2(n + l)K n+1 {H) 
m 

+ ((2n + 1 W - ^0f)K n (H) + \AnK n _i(H)] (4) 
L n (H) = -^^K n (H) - h{n + l)K n+1 (H) . 

If h is set to zero in Eq. (4) and H replaced by a numerical value E, one recovers the 
difference equations for the corresponding classical equation mentioned above. Our starting 
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Eq. (2) can be compared with the operator equation [F(x,p), H] = ih, a quantum version 
of {4>/uj(H), H} = 1, solved by C. Bender [Q with a different type of ansatz (which requires 
negative powers of p or x). 

We now study the matrix element < >= (hKd/dx + L) < x\E >, where L and K 

are short notations for Y^=o L n (E)x 2n+1 and Yl^Lo E- n (E)x 2n respectively. In the following, 
the dependence of L, K and A on E is implicit and primes denote ^-derivatives. Eq. (2) is 
satisfied for an arbitrary potential provided that 

JL K » + ^ + KA = (5) 
2m m 

h 

L" - 2K'(E -V) + KV' + LA = . 

In the case of the potential given in Eq. (1), these equations can be solved using the 
linear recursion relations given by by Eq. (4) with H replaced by E. As a result, we get 

< > as a function of E, A(E), K and K\. If both K and K\ are zero, then 

< >= 0, consequently we require at least one of them to be non-zero. In order to fix 
the ambiguity mentioned above, we shall impose that the first non-zero K n to be 1. One 

2 

can prove rigorously that \K n \ < Ci(C2) n (n\)~3 with C\ and Ci independent of n. This 
bound implies that L and K are entire functions of the position seen as a complex variable. 
This allows us to controllably approximate these functions in terms of their truncated power 
series. Note that in the classical case, the factorial suppression is absent and the functions 
L and K have a common finite radius of convergence which reflects the existence of turning 
points. 

We are now in position to construct formal solutions of the time-independent Schrodinger 
equation. First we notice that 

L 2 - h(L'K - K'L) - 2mK 2 (V - E) = K (h 2 Ki + hmAK + 2mK E) . (6) 

Recalling Eq. (4), the equality clearly holds at x = 0. In addition, Eq. (5) implies that the 
derivative of the LHS is zero, consequently it holds for any x. Adjusting the constants K 
and Ki in such a way that the LHS of Eq. (6) is zero, and dividing by K 2 (temporarily 
assuming that K ^ 0), we obtain that L/K satisfies the Ricatti equation 

h(^)' = (|) 2 - 2m(V - E) . (7) 

A detailed study of Eq. (5) shows that Eq. (7) is also satisfied near a zero of K. Parametriz- 
ing the wavefunction as 

< x\E > oc e~^^" dv ^) , (8) 

we find that Eq. (7) is the Ricatti form of the time independent Schrodinger equation. Note 
that Eq. (8) implies that A\E >= 0. We still need to specify the conditions under which the 
RHS of Eq. (6) is zero. After fixing the arbitrariness in the coefficients as discussed above, 
there remains two possibilities. The first one is K = 1 and K x = — 2rr ^° (4p + E) which 
corresponds to an even eigenfunction < x\E >. The second possibility is Kq = and K\ — 1, 
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which corresponds to an odd eigenfunction, and for which L/K has automatically the —H/x 
singularity at the origin. In both cases, Eq. (4) define uniquely L(x) and K(x) given E and 
A. From the uniqueness of the solution of the Ricatti equation, given a condition at x = 
(L/K = in the even case and K/L = in the odd case), L/K is indeed A- independent as 
one can check order by order in the expansion of L/K near x = 0. A particularly convenient 
choice is A = 0, because in this case Eq. (5) and (8) imply that | < x\E > | 2 oc |X(a;)|. 
We shall now use the formal solution of the Schrodinger equation to find sharp upper and 
lower bounds on the energy levels. 

Our basic tool to find accurate upper bounds on the energy levels will be the theorem 
proven in Ref. [5] that for a Sturm-Liouville problem the n-th eigenfunction divides the 
fundamental domain into n parts by means of its nodal points. For the problem discussed 
here, the zeros of < x\E > are the poles of L/K which is seen most easily by picking 
A = 0. Consequently, if K has more than n zeros at finite x (which are not zeros of L), 
then E > E n . Furthermore, if E is decreased continuously, the largest zero of K increases 
continuously. When E reaches E n , a pair of zeros disappears at infinity. One can then 
monitor the "entrance" of the zeros in the region x < a while E increases, by finding E 
such that K(a, E) = 0. This can be done using Newton's method with a an appropriate 
truncation in the expansion of K. One can then check if the existence of the zero can be 
established despite the errors due to the truncation (which can be estimated using the bound 
mentioned above). 

Lower bounds can be found from the requirement that the wave function < x\E > is 
normalizable. Due to the fact that < x\E > has a definite parity, we will restrict the 
discussion to positive x part of the (x,L/K) plane. This half plane can be divided into a 
region where L/K increases and a region where L/K decreases. The boundary between these 
two regions is characterized by (4)' = which by Eq. (7) implies (4) = ±^/2m(V — E). 
For this reason we call this curve the 'WKB curve". For the potential of Eq. (1), one finds 
that if a trajectory crosses the WKB curve L/K continues to decrease when x increases. In 
other words, the WKB curve is the boundary of a sink. Using bounds on (L/K)' coming 
from Eq. (7) we can prove that whenever the trajectory crosses the WKB curve, < x\E > 
defined by Eq. (8) is not normalizable. Consequently, if < x\E > has n nodes and if the 
corresponding L/K ultimately decreases (when x becomes large enough), then E < E n . 
Note that Eq. (7) shows that the poles of L/K can only be simple and have a residue 
—h. This implies that the Dunham condition [[F] used in semi classical calculations || is 
automatically satisfied. This also implies that for the potential of Eq. (1), a normalizable 
wave function cannot have a zero in the classically forbidden region (because L/K could 
never reach the positive part of the WKB curve for x larger than the location of the pole). 

In summary, when E is sufficiently close to E n and x sufficiently large, L/K follows 
closely the trajectory of the positive part of the WKB curve. When a certain value of x is 
reached L/K depends sensitively on small changes in E. A small increase in E creates an 
additional zero of the wavefunction, a small decrease forces L/K to cross the positive part 
of the WKB curve and to reach the negative part of it. This allows us to find sharp bounds 
on the energy levels. The only problem which remains is the control of the round-off errors. 
For a usual double precision computation, this is a serious issue, however since the linear 
recursion formula of Eq. (4) requires a number of operations which only grows linearly 
with the maximal order calculated, we can use "slow" computational methods involving 
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a very high precision. This can be implemented, for instance with Mathematica, using 
the instruction "SetPrecision[...,100]" for numbers set with a precision of 100 digits. This 
method has allowed us to obtain the wave function with very good precision, at large x, far 
beyond the classical turning point, i.e., where the lowest order WKB approximation works 
well. 

Proceeding this way, we have calculated the first 10 energy levels in the case m = 1/2, 
ujq = 2 and A = 1/10. The results are displayed in Table 1 with 30 significant digits. 
These numbers have been obtained by keeping 400 terms in the expansion of L and K 
and restricting the calculation to the interval \x\ < 7.5 The starting precision was 100 
digits. The difference between the upper and the lower bounds were required to be less 
than 10~ 32 . These calculations have been performed independently using Mathematica and 
Maple. Our numerical results are in agreement with the existing literature summarized in 
Ref. [6] and where numbers up to 15 significant digits can be found. Note that the numerical 
precision on the lower bounds obtained with Newton's methods (which can estimated from 
the Mathematica command "Precisionf...]") decreases approximately linearly when the level 
increases as shown in the last column of Table 1. It is clear that these round-off errors are 
much smaller than the theoretical precision achieved. More generally, a preliminary analysis 
indicates that the enterprise of calculating a very large number of levels with a very large 
precision does not face prohibitive (i.e, exponential) growth of computer time. If this is 
effectively the case, we could say that the quantum anharmonic oscillator is "numerically 
integrable" . 

A more satisfactory outcome would be to find an operator which would be the analog of 
the classical action which satisfies the Poisson bracket relation {1(H), e 1 ^} = —ie 1 ^. At the 
quantum level the corresponding relation implies an equally spaced spectrum. Equivalently, 
if we had an analytical expression for A(E) corresponding to A being the minimal creation 
operator (where E is replaced by E n and hA(E) by E n+ % — E n in the matrix elements), we 
could calculate recursively the energy spectrum. In both cases, it would mean that we would 
have at hand an implicit closed form expression for the energy spectrum. Despite interesting 
attempts, || such an expression has not been found and discovering it is a challenge for the 
future. 

It is a pleasure to thank T. Allen, A. Bhattacharjee, C. Bender, W. Klink, Y. Nambu, 
W. Polyzou, G. Payne, V. Rodgers, J. Schweitzer, D. Speiser and J. Weyers for valuable 
conversations and comments. 
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Table 1 



Tl 




s.d. 


o 


1 06528550954371768885709162879 


95 


1 


3 30687201315291350712812168469 


93 


2 


5 74795926883356330473350311848 


89 


3 


8 35267782578575471215525773464 


87 


4 


11.0985956226330430110864587493 


84 


5 


13.9699261977427993009734339568 


81 


6 


16.9547946861441513376926165088 


79 


7 


20.0438636041884612336414211074 


77 


8 


23.2295521799392890706470874343 


74 


9 


26.5055547525366174174695030067 


72 



Table Caption 

Table 1: The first ten energy levels E n and the number of (numerically) significant digits of 
the upper bound on E n (s.d) in the case A = 1/10, m = 1/2 and uj = 2. 
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